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Abstract 

We consider a stochastic model of clock synchronization in a wireless 
network consisting of N sensors interacting with one dedicated accurate time 
server. For large N we find an estimate of the final time sychronization error 
for global and relative synchronization. Main results concern a behavior of 
the network on different time scales tjy — > oo, N — > oo. We discuss existence 
of phase transitions and find exact time scales on which an effective clock 
synchronization of the system takes place. 

Keywords: Clock synchronization, time to synchronization, timeliness in real-time 
systems, wireless sensor networks, multi-dimensional Markov process, self-organization, 
phase transitions 

1 Introduction 

For many years distributed systems are a constant source of challenging research 
problems. Clock synchronization is one of the most well-known and widely dis- 
cussed subjects. In distributed systems there is no global clock. Clocks of different 
components tick at different rates. But to achieve an efficient parallelization in 
their common job all components of the distributed system need a common no- 
tion of time. This is very important for real-time systems where scheduling and 
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timeliness play a crucial role (see, for example, [7J [32]). Real-time applications 
must guarantee response within strict time constraints. Similarly, Wireless Sensor 
Networks (WSNs), having their own specific character among distributed systems, 
usually meet many real-time requirements. In particular, energy saving algorithms 
in WSNs would be impossible without a global consistent timescale |27] . Nowa- 
days complex real-time embedded systems are distributed and have many active 
components. Modern WSNs also consist of very large number of sensors (|8| 124]). 
Hence an important requirement for clock synchronization in distributed systems 
is the scalability of corresponding algorithms. 

In this paper we consider a mathematical model of a large WSN communicating 
with one dedicated server of accurate time. The model is based on a multi- 
dimensional stochastic process. The network consists of N sensors equipped with 
non-perfect clocks. By using message timestamps sensors share information about 
their local times and receive accurate time from the dedicated server. It is assumed 
that all messages are sent at random time moments. Main focus of the paper is the 
asymptotical behavior of the network in the limit when the number N of sensors 
tends to infinity. 

This article is organized as follows. The model is defined in Section |2j It de- 
pends on a few parameters. The network is assumed to be symmetric. Here we 
do not touch many problems (synchronization protocols, energy saving optimiza- 
tion etc.) which are very important for practical use of WSNs. Such questions 
were widely discussed by many authors [30| [9j [271 El EH]. The present article is 
a mathematical paper. We intentionally make our model as transparent as pos- 
sible to keep this paper short and to illuminate main results about existence of 
different phases in the evolution of the network (Theorems 1-3 in Section [3]). It 
appears that the proposed model has almost explicit solution: in all asymptotical 
expressions we can play with parameters and discover many interesting phenom- 
ena related with large networks. In Section [6] we discuss possible generalization 
of the present results to more general models. We do not seek for direct practical 
implementation of our results. But we believe that such results are interesting not 
only from theoretical viewpoint. 

Our probablistic technique is similar to that was recently used to study math- 
ematical models of multi-processor parallel computing |2T1 [H HH [13] and col- 
lective behavior in abstract particle systems with synchronization-like interac- 
tion [121 Cni EE EJ [161 HH [19] . Our stochastic models can be considered also as 
special classes of self-organizing systems [2H1 El 1231 G5] - 
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2 Mathematical model 

There is a wireless network with N + 1 nodes. Each node i has its own clock. 
Let Xi(t) be the value of clock % (the local time at node %). The physical time is 
denoted by t G R+. Nodes 2, 3, . . . , iV + 1 correspond to sensors. We assume that 
clocks of the sensors are identical but not perfect and progress in the following 
way: 

x j (t) = x j {0)+vt + aB j (t), j = 2, . . . , N + 1 . (1) 

Node 1 is a time server, we assume that its clock is perfect. Clock 1 reports a local 
time X\(t) which is a linear function of the physical time: 

xi(t)=Xi(0)+rt. (2) 

The constants r and v are positive, they are called frequencies of the corresponding 
clocks (see, for example, [30] ) - In general, these frequencies are not equal: r ^ v. 
The parameter a > in (pQ) corresponds to the strength of a random noise related 
with imperfect clock of a sensor. For simplicity the random noise in flTJ is modelled 
with independent standard Brownian motions 

(Bj(t), t>0), j = 2,...,N+l. 

This white noise assumption is usual for many oscillator clock models [32| [23]. 

As it was mentioned in Section [T] (see also [28J) the problem of clock synchro- 
nization is critical for the proper work of WSNs. To synchronize local times of 
a pair of sensors the most natural solution is to send time-stamped messages be- 
tween them. After reading a new coming message the receiver ajusts its clock to 
the local time of the sender recorded in this message. Once being synchronized 
the clocks of the pair will diverge immediately due to assumptions made in flTJ 
and Q2). 

We continue with the formal definition of the model. We shall call ([I])-® a 
free dynamics. The free dynamics generates independent evolutions at nodes of 
the network. Now we add some special interaction between nodes. Namely, at 
random time moments each node sends messages to other nodes. Our assumptions 
are 

• with the rate a > the server node 1 generates a message containing infor- 
mation about the current value of x\ and sends this message to one of the 

sensors which is chosen randomly with probability — ; 

• independently each sensor, with the rate > 0, generates a message about 

its local time and sends it to another sensor which is chosen randomly with 
1 

probability 



N-l 
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• messages reach their destinations instantly (there are no transmission de- 
lays); 

• if sensor j receives a message from sensor % at some (random) time r then 
clock j is immediately ajusted to the value of clock i: xj(t + 0) = Xi(r); 

• between receiving of subsequent messages sensor nodes evolve according to 
the free dynamics ([T]) ; 

• server node 1 always follows the free dynamics (j2J). 

Let us remind some terminology used above: a sequence of events = r < • • • < 
r n < ■ ■ ■ generated with rate 5 > is called also a Poisson flow of intensity 5. 
It means that intervals between events r n+ i — r n are independent exponentially 
distributed random variables with mean 7 _1 : 

P (t„+i - r n > s) = e~ 5s , E (r n+1 - r n ) = 6' 1 . (3) 

Hence we defined the multi-dimensional stochastic process 

x (t) = (x 1 (t),x 2 {t), . . .,x N+1 (t)) , t > 0, 

which appears to be a continuous time Markov process with values in M 7V+1 . This 
process is an idealized mathematical model of the homogeneous wireless sensor 
network with dedicated accurate time server. 

It should be noted that non-Markovian models of WSNs can be also consid- 
ered in the framework of the present paper. We postpone discussion on possible 
generalizations to Subsection 14.21 

To estimate desynchronization in the network it is convenient to consider the 
following functions on the configuration space M Ar+1 : 



1 N+l 

D<kX) := (N-DN foa ~ x n) 2 - (5) 

The first function corresponds to deviations of sensor's clocks from the accurate 
time X\. The second function is related only to internal inconsistency between 
the sensor's local times x 2 , ■ ■ ■ , Xn+i- This function D(x) may be more useful 
for causuality-based real-time models where the right ordering of events is more 
significant than the global synchronization to physical time (see, for example, [31]). 
Since the both functions R(x) and D(x) are averagings of {xf. — Xj) 2 , the squares 
of offsets between clocks, the true sense of a time synchronization error have their 
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square roots a/ R(x) and yjD{x). The physical dimension of ■sj R(x) and D{x) 
is time, i.e., they are measured in seconds. 

Note that for any t values R(x(t)) and D(x(t)) are random. It is more conve- 
nient to deal with their expectations: 

R N (t) := ER(x(t)) and D N (t) := ED(x(t)) . 

3 Main results 

Let x(t) be the Markov process introduced in Section [21 

In Subsections I3.1H3.3I we present various results on asymptotical behavior 
of the functions -Rjv(t) and -Djv(t). These results will provide us with detailed 
information on collective behavior of the network in the limit when N — > oo. Here 
we give one particular corollary of our theorems which is a good illustration to an 
approach followed in this paper. 

Corollary. Assume that u/r. For simplicity take Xi(0) = ■ ■ ■ = xtv(O) = as 
initial state of the network. Let s > and 7 > be parameters of a new time 
scale. Namely, we put t = siV 7 and look after the function D^{t) when N — > 00. 
It appears that 

Djv(siV 7 ) ~ C(s,j) 

where 

(7, 7<i 
0(7) = < 3 7 - 1, \ < 7 < 1, 
[2, 7 >1. 

The function C(s, 7) = C(s, 7; <r, a, /3) > depends on parameters a, a, (3 and is 
increasing in s for fixed 7. 

The choice £ = sN 1 means that we consider a new time unit iV 7 and the "size" 
of configuration X2, ■ ■ ■ , x^ on this scale is of order N^ 1 ^ 2 . The conclusion is that 
on different time scales the large network shows very different types of behavior. 
Note that the function ^(7) is not smooth. This situation is very like the phase 
transitions phenomena in models of statistical physics. We give more details on 
time scales in Subsection 14.61 

To understand in what sense the algorithm of Section [2] drives the network to 
synchronization one should answer several questions. 
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3.1 N is fixed, t — > oo 



The first question is a long-time behavior of the stochastic process 



x{t) = (xi(t), . . .,x N+1 ) . 



Since r > and v > it is quite clear from definition of the model that Xj(t) — > oo 
(t — > oo) for each node j. To have more detailed information on the process we 
put a moving observer at the point Xi(t). From the viewpoint of this observer 
states of the sensor nodes 2, . . ., N + 1 are given by the vector 



It follows from the general theory of Markov processes (see Doeblin condition 
in [3]) that y(t) has a limiting distribution as t — > oo. Similar situations were 
discussed in [T5j [TS] for different synchronization models. Hence our first conclu- 
sion is: for fixed parameters a, a > 0, (5 > and iV the Markov process y(t) 
is ergodic and converges to its equlibrium as t — > oo. This limiting distribution 
(which is a probability measure on R ) can not be obtained in an explicit form. 
Nevertheless, some important functionals of the limiting distribution can be found 
explicitly. Note that R(x) = R{y) and D(x) = D(y) where x and y are related 
by (JSD- Thus the clocks Xi(t) of sensor nodes are synchronized in the following 
sense: mean values of time synchronization errors stabilize and do not change in 
time anymore. Namely, assuming that N is fixed and t — > oo we can prove the 
following statement. 

Theorem 1. We distinguish two cases: the skew v — r of a sensor's clock relative 
to the accurate time is zero or nonzero. 

Case 1: Assume that v = r. Then 



yit) = (yz(t), yN+i(t)) , yj = Xj - xi- 



(6) 



R N {t) -> fl^(oo) : 




a 



D N (t) -> D%\oo) 



2 a 2 

a + (3 



N . 



Case 2: Assume that v ^ r. Then 



R N {t) -)• ^ 2 (oo) : 




N 2 . 



D N (t) Dtf(oo) 




N 2 . 



a (a + 0) 
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As usual we write here and below fjy ~ gN iff li m — = 1- 

N 9n 

Now we come to a natural and very important question: how many time we 
should wait until the networks of "size" N will be synchronized ? We devote to 
this question the next Subsections I3.2H3.41 

Comparing Case 1 and Case 2 in the above theorem we see that result of the 
global synchronization of sensor's clocks is worse in the biased case v ^ r. Indeed, 
condition v ^ r means that there is a systematic error in client's clocks which 

appears to be more essential (a resulting "clock offset" is of order \J R^ 2 (oo) = 
O(N) ) than in situation of a pure random noise errors under condition v = r 

(when the "clock offset" is of order ^ R^(oo) = C^iV 1 / 2 ) ). 

Note also that in Case 2 the limiting values i?^ 2 (oo) and D^ 2 {oo) do not 
depend on er, the noise parameter of the sensor's clocks. Moreover, these values 
give the right asymptotics even for a = 0. This means that the so called locked 
synchronization [28J is not possible in our model due to the stochastic nature of 
the message-passing algorithm. 



3.2 N — > oc, £ — > oo. Phases of synchronization 

Synchronization in networks with large (or growing) number of nodes is of spe- 
cial interest. Our approach is to consider limits when both the number of client 
nodes iV and the time t grow to infinity. More precisely, we shall consider se- 
quences (AT, tjv); where physical time t = t/v is some increasing function of N. We 
shall see below that for different choices of tpj network exhibits different asymp- 
totical behavior. One can say that a large network passes different phases on its 
road to synchronization. For our model there exist several time scales t = tff of 
qualitatively different behavior. We always assume that iV — > oo. 

To begin we assume in the present subsection that initial distribution of the 
local clocks is such that sequences {-Rat(O)} and {Dn(0)} are bounded in iV: 

sup-Rjv(O) < +oo, sup D N (0) < +oo. (7) 

N N 

In Subsec. 13.31 we shall discuss the situation when {i?7v(0)} and {Dn(0)} are 
unbounded. As in Subsect. 13. II we present results separately for Case 1 and Case 2. 

Theorem 2. Let assumption hold. 

Case 1: v = r — zero skew. 

PI. tjy/N — > (phase of initial desynchronization) : 

( R N (t N ) \ (a 2 \ 
V D N (t N ) ) ~ (, 2a 2 ) tN 
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P2. tx/N — > c > (phase of effective synchronization) : 

where Id is the identity map in M 2 , 

M ■— ( ~ a ° 
• V 2a -2(a + 0) 

P3. tw/N— 7- +oo (phase of final stabilization): 

( R N (t N ) \ ( a 2 /a \ / flQ^oo) 

V D*(fjr) ; ~ ^ 2 ( r 2 /(a + /3) J ~ ^ D^(oo) 

Remark 1. By using definition (j8J of the matrix M it is easy to check that M 
has two distinct negative eigenvalues: Ai = —a, A 2 = —2(a + /3). Moreover, 
(-cM)- 1 (Id - e cM ) ->■ Id as c +0, e cM and 



-M)- 



l\ ( 1/a 
2 J \2/(a + /3) 



From this remark we see that phase P2 joins smoothly asymptotics of phases PI 
and P3. 

Remark 2. One can easily calculate R^(t) and D^f(t) in the simplest case v = r, 
a = P = (no bias, no synchronizing interaction). It appears that -Rjv(i) = 
a 2 t and -Djv(£) = 2cr 2 t. This observation brings us to the following explanation 
of the phase 1: on time intervals [0, o(N)] the influence of synchronizing jumps 
(clock adjustments) is negligible with respect to the impact of the random noise 
of the free dynamics. On the time scale of phase P2 (t^ = cN) the cumulative 
effect of individual synchronizing adjustments become of the same order as the 
effect of random noise, so we obtain the non-trivial dependence of Rn^n) and 
Dn^n) on c. On phase P3 the values of Rn^n) and Dm^n) correspond to the 
synchronized network. It means that the effective synchronization takes place on 
the times of order N. 

Theorem 3. Assume that condition ffi) holds. 
Case 2: v 7^ r — nonzero skew. 
PI. tx/N — > (phase of initial desynchronization) : 

RnM ~\(v-r) 2 t 2 N . 
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Pla. // 



Plb. // 



1~N 



N 



— > 0, then D N (t N ) ~ 2a 2 t N . 



<N 
Dn^n 



— > ci, c\ > 0, i/ien 
1 



2a 2 + -a(v - r) 2 c\ ) t N . 



Pic. // -^L^oo but ^^0, tfien ~ \a{v - rf tl, / N . 



P2. t^/N — i- c > (phase of effective synchronization): there exist functions 
h R and ho not depending on N such that 



Rn^n) 
D N (t N ) 



h R {c) 
h D (c) 



[v — r) t 



2 + 2 

iV 



P3. ttf/N —7- +oo (phase of final stabilization): 



( (v — r) 2 \ 



Rn^n) 
D N (t N ) 



Or 



[y — r) 2 
V a(a + /3) ) 



N 2 




Remark 3. Explicit form of the function h R and hn is given in (136|) . see Section H75J 
It appears that 



h R (c) 
h D (c) 



1/2 

ac/3 



(c^+0), 



(c +OOJ. 



h R (c)c 2 \ a- 2 \ 

h D (c)c 2 J \ a~\a + /3)- 1 J [C 

Hence the phase P2 "continuously" joins asymptotics of phases Pic and P3. 



3.3 Situation of a "big initial disorder" 

As in preceeding subsection we consider here different time scales t = t^ assuming 
that t]y — }■ ooasiV — y oo. Our goal is to discuss what happens with sychronization 
phases if the sequences R^(0) and D^(0) grow to infinity as N — > oo. This 
problem is not too complicated and can be solved explicitely for any concrete 
assumption about initial disorder, i.e., about increasing rate of Rn{Q) and Dn{0)- 
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In this paper we shall not give universal and general answer to this question since it 
cannot be presented in short and transparent expressions. To avoid cumbersome 
formulae we just describe the answer in a few general words. There exist a so 
called initial disorder decay interval (0,t° N ). The function t° N depends on R^(0), 
Dn(0) and all parameters of the model. Under above assumptions t° N — > oo as 
N — > oo. If tjv — °(^Ar) 5 then the values Rn{0) and Dn(0) enter in asymptotics of 
Rn^n) and D^(t^), but they don't enter in these asymptotics for such that 

Given a function t° N , one should compare orders of t° N , N 1 ^ 2 and N to see 
intersections of the initial disorder decay phase (IDDP) with the phases PI, P2 
and P3 of Subsection 13.21 The IDDP dominates over any of the phases PI (Pla- 
Plc), P2 and P3. So depending on the initial disorder in the network one can see 
the following different sequences of consequtive phases in evolution of the network: 
IDDP-P1-P2-P3, IDDP-Plb-Plc-P2-P3, IDDP-Plc-P2-P3, IDDP-P2-P3 or 
IDDP-P3. 

3.4 Collective displacement from the etalon time 

In Subsections I3.lf l 3.3l we studied R(x) and D(x) which are quadratic functions 
of the clock configuration x G M^" 1 " 1 . For completeness we consider here a linear 
function d(x) which gives some information on a collective displacement of the 
sensors clocks X2, ■ ■ ■ , %n from the clock of the server x±. Define functions 



d N (t) = Ed(x(t)) , d N : R + -> R. 
Theorem 4. Assume that sup Ar cZ A r(0) < +oo. 

Case 1: v = r — zero skew. It appears that djsf(t) — >• as t — >■ oo for any fixed N. 
Case 2: v ^ r — nonzero skew. For fixed N 




JV+l 




(9) 



and 



d^it) {y — r)a N as t — > oo. 



For time scales t^ — > oo (N — > oo) the following statements hold: 
o ift^/N — » 0, then djv(^Ar) ~ (v — r)tw , 
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o if t N /N — >■ c> 0, then 

rfjv(^jv) ~ (1 — exp (— etc)) (u — r)a _1 iV ~ - 
o ift N /N — > +oo ; then d N (t N ) ~ (v — r)a~ l N . 

The proof of this theorem is similar to proofs of Theorems [TH2] (see Section HJ) 
but technically it is much easier. So we omit it. 

Remark 4. Let us consider the degenerated model with a = 0. It is easy to prove 
that 

d N (t N ) - d N (0) = (v- r)t N . (10) 

Comparing this with Theorem H] we see that the formal limit a — > does not turn 
results of Theorem H] into (ITD]) . Hence we conclude that limits iV — > oo and a — )■ 
do not commute. 



— exp(— ac) 

ac 



r)t 



N 



4 Proofs 

4.1 Conditional averaging 

Let = t < ■ ■ • < r n < ■ ■ ■ be sequence of time moments when messages are sent 
(received). It follows from definition of the model (Section[2]) that {r m — T m _ 1 }^ =0 
are independent i.d. r.v. having exponential distribution with mean (a + iV/3)" 1 . 
In the sequel we will refer to this condition on {r n } as to a Markovian assump- 
tion. 

Let U t = max{m : r m < t}. To get R^t) and D^(t) we will calculate the 
chain of conditional expectations as follows 

E(-) = E (E (E (-Kr,}^) |n,)). 

Let / = f(x) be some function on the configuration space R 7V+1 . Introduce 
notation 

/ (n) = E (f(x(r n )) | fa}^) , n = l,2,.... (11) 

Hence f^ is a random variable functionally depending on the sequence {t?}^. 
In other words, to any function / = f(x) we put in correspondence / i— > {f^} 
the sequence of random variables {f^}. 
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Denote 



V{x) = ^ ) X [ , V: R N+1 M 2 , 



D(a?) 

where functions R(x) and -D(a;) are defined in (SD — dSD • Recall definition of the 

function d(x) from (jUJ). In the next subsection we study random sequences 

and {c/( n )} corresponding to the above functions V(x) and d(x) according to ( ITT]) . 

4.2 Recurrent equations 

r , / \ f s 

In this subsection we obtain difference equations for V w = I 1 and d^ n > . 

To present this result we need some convenient notation which will be used till 
the end of the paper. Namely, we introduce the following reals 



a (3 
v-r, a N = — , (3 N 



N 1 ' N-l' 

u = 2b, 8n = a + N(3, k N = 1 — 5^ (%n 

and two-dimensional vectors 

<li = ° 2 ( o ) ' ^2 = (v ~ rf ( J , q = w ' 



:i2i 



' HU 



We consider also a 2 x 2-matrix 



2« w -2(a N + /3 N ) 

and define a linear operator = Id + S^L^, where Id is the identity map in IR 2 . 
Lemma 1. The following recurrent equations for V^ n ' and d^ hold 
rfW = k N {d {n ~ l) + 6A n ) 

V™ = K N + A 2 ? 2 + A nQl + A n d^q ) 

where A n = r n - r n _i . 

Proof of Lemma [1], The dynamics of the process x(t) consists of two parts: 
free motion and pairwise interaction between nodes of the network. Namely, the 
interaction is possible only at random time moments < T\ < T2 < ■ ■ ■ and has 
the form of synchronizing jumps: at time r n a pair of nodes 

i = l,2,...,N, j = 2,...,N, i^j, 
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is randomly chosen, and the value of clock j jumps to the value of clock i: 

(xi,Xj) -)> (xi,Xi) . (14) 

Inside intervals (Tk,Tk+i) components of the process x(t) move according to the 
free evolutions (pQ) and (J2j) driven by independent Brownian motions. Recall that 
the pair mentioned above corresponds to a time-stamped message sent by 

clock % to clock j at time r n . It is easy to see that any pair is chosen with 

probability a^d^ 1 and any pair i > 1, is chosen with probability fijqo'jj 1 - 

Keeping in mind (1141) we introduce a family of maps Suj\ : ~R N+1 — > ~R N+1 , 

'■ { x li • • • ) x N+l) ^ (^1) • ■ • j x N+l) J 

where a^- = Xi, x' k = Xk, k ^ j, and define a map- valued random variable S 
such that 

P{S = S {lj) } = h 5 \> < = 1 ' (15) 

So Lemma [T] will immediatelly follow from the next two lemmas. 

Lemma 2. For any x G IR Ar+1 

Ey(Sx) = K N V(x), Ed(Sx) = k N d(x). 

Denote an auxiliary stochastic process z — z(t) evolving according to the free 
dynamics flU)-©- 

Lemma 3. For s <t 



E (d(z(t))\z(s)) = d(z(s)) + bA 

E (V(z(t)) | z(s)) = V{z{s)) + A 2 • q 2 + A • ft + A • d(^)) ?0 
where A = t — s. 

The proofs of Lemmas |2] and [3] are quite straightforward and very similar to 
proofs of corresponding lemmas in [T6| [T8| [19] . They are omitted here. 

It should be noted that, in fact, Lemmas [2] and [3] are valid under weaker condi- 
tions than the Markovian assumption. Namely, these lemma are true also for the 
following semi- Markovian assumption: messages are sent at epochs {t„}^ 
where intervals {r m — T m „ 1 }^' =0 are independent random variables with identi- 
cal continuous distribution; at time moment r n the pair (i,j)= (sender, receiver) 
is chosen independently according to the map-valued random variable S defined 
in (US}. 

Briefly speaking, under Markovian assumption the sequence {t ti }°° =0 is a Pois- 
son flow but under semi-Markovian assumption {r n }^ =0 is a general renewal pro- 
cess [21 S]. 
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4.3 Decomposition into sum over diagrams 

Having equations of Lemma [T] we can use them to express via V( n -V and 
so on. Proceeding recursively we get a sum with large number of summands. 
We want to organize these summands in some proper way to be able to evaluate 
as an explicit function of V^. To do this we use a general approach known 
as decomposition into sum over diagrams. We start from description of a set of 
diagrams corresponding to our specific task. 

Admissible diagrams. Fix some n. Let us define a set of admissible diagrams 
Q(n) of order n. We say that an oriented graph G belongs to the set Q(n) iff G is 
a path 

G = (v , t>i, . . . , v r ) , 1 < r < n + 1, (16) 

with vertices v , Vi,...,v r G Ai n and oriented edges v^i —> v iy i = l,...,r, 
satisfying to the following conditions. 

a) Vertices of G are labeled by pairs (c, d) belonging to 

M n = { v = (c, d) : (c, d) e 4 n) xS D U {("t", 2)} } 

where c G := {0, 1, . . . , n — 1, n} and d G So '■= {0, 1, 2}. 

b) The initial vertex Vq of any path G is ("t",2). For % > 1 the vertex t), has 
the form (n + 1 — i, Zj G Sx>. The final vertex may be 

{(n + 1 — r, 0), if r < n, 
(0, d), if r = n + 1 where ci G S*/). 

c) The first edge vq — > v\ of G has the form ("t",2) — >■ (n,V), I' G S 1 /), and 
other edges (t>j — >■ f 1 < i < r — 1, have the form (n + l — i, k) — > (n — i, k + i) 
where < l^. There are no edges of the form (m, 0) — > (m — 1, 0). 

For the path (fT5|) we use notation 

V G = {^o, v 1 , ...,v r } , £ G = {v^i -)> | i = 1, . . . , r} 

for the set of vertices and the set of edges of G. For reasons which will be clear 
below, sometimes we will remap elements of the set Sd as follows 

0< — ► "1", U — > "d", 2< — > "V". (17) 
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Contribution of a path. Let f — (tl, . . . ,r„) G M + and t > are such that 
< 7~i < ■ ■ • < r n < t. We define Jq = Jg(t, t) 6 K 2 , a contribution of the path G 
(see (TTSj) ). as follows 

\a£V G be£ G / ordered along the path 

where functions p a = p a (r,t) and 717, = 7^(7, t) are contributions of vertices and 
edges defined below. Namely, we put 



Pv — 


1 










if f = 


(m,2), 




ry(o) ) 


if v 


= (0,2) 


: H 




if v = 


(to, 1), 


Pv = l 


If' 


if v 


= (0,1) 




I 1 ' 


if f = 


(m,0), 






if v 


= (0,0) 






m — 1, 


...,n. 


and 
















Kb = fi,i'( 






m — 1, . 


..n+l, 


Tn+1 = * 






for the edge b - 


= (m,Z) 


(m 


-i,0- 


The functions /^/ are 


defined 


as follows: 



f 22(3) = Id, f w (s) = s 2 q 2 + sq l5 / u (s) = 1, 

/21(a) = sqo> / 10 (s) = sb. 

Proposition 1. The following expansion holds 



E (v(*(f)| fa}^) = £ J G , 

G€6(n t ) 

where Q{n) is the set of admissible diagrams of order n (defined above) and Jq e M 2 
is t/ie contribution of a diagram G. 

The proof of Proposition [TJ is just a careful development of the recurrent equa- 
tions of Lemma [TJ The only thing one should pay attention is the last inter- 
val (rn t ,t] of the total time segment [0,t]. There is no synchronization jump at 
time t. Hence we need to apply for this interval only Lemma [3j This explains 
assignment p VQ = 1 for the root vertex vq = ("i", 2). 
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4.4 Evaluation of the functions 

Below we use agreement (II 7p . For any path G G Q(n) denote by n\ = n\{G) 
the number of vertices of the form (m, "d"), m < n, and by n 2 = n 2 {G) the 
number of vertices of the form (m, "V"), m < n. In is clear that nx,n 2 > and 
^1 + n 2 < n + 1. Now we decompose the set of admissible diagrams Q(n) into the 
following nonintersecting subsets 

G(n) = Q (n) U g 2 (n) U g 10 (n) U £ n (n), 

g (n) = {Ge g{n) : m = 0, n 2 = n + 1} = {G } 
2 (n) = {GE d{n) : m = 0, < n 2 < n} 
gio( n ) = {G E g(n) : 1 < ni < n + 1, n 2 = n + 1 — %} 
£n(rc) = {G G <?(n) : 1 < ni < n, n 2 > 0, n\ + n 2 < n + 1} 
The subset go(n) consists of a single path 

Go = (("*", "V"), {n,«V), (n-l,«V"),..., (0,"V)) 
and contribution of this path is 



JG (r,t) = Y[ n JJ Pa 



K n N V {0) . 



G=G 

For G G g 2 (n) the final vertex u r is (n — n 2 ,"l") and 

J G (f,t) = K^A 2 n+1 _ n2 q 2 + K%A n+l _ n2qi , 
where A k = r k — Tk-i- If G G ^io(^) then v r = (0,"<i"), r = n + 1, and 
J o (7' l t) = J K5+ 1 - B1 ^- 1 A Bl dP , )q 
For G G ^n(n) the final vertex v r is (n — n\ — n 2 ,"l") and 

Jg\Ji^) = Kjy k 1 ^ A n+ i_ n2 A n+ i_ ni _ n2 bq 
So we have 



E (V(x(t) | {75}°!,) = f 



E + E + E + E I*. 

Ge£?o(n t ) Geg 2 (nt) GeSia(n t ) Geen(n t ) 
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taking conditional expectation E (• | ILJ of the both sides of this decomposition we 
obtain 



E(V(x(t)\U t = n) = (E+E + E + E| J < 

.Geg (n) Gag 2 {m) Geg 10 (n) Geg n (n) 



G 



n 

K%V® + (n 2 4 2 il-n 2 q 2 + ^Sl-n^l) + 

n 2 =0 



vn 2 =0 
' n—1 n— n 2 



4- V V /r n2 A- ni s (1,1) bn 
> \ / j / j N N n+1— n 2 ,n+l— ni— n 2 / "4' 

^n 2 =0 m=l 







where 

4 m) = si m) (t,n) = E((A fc r|II t = n) > m = l,2, 
s 5 1} = *S5 1) (*,n) = E(A i A i |II t = n) l 

It is easy to check that in general semi-Markov case (Subsection 14. 2ft 

s hfi ' n ) = s i 2 ',h (* > n ) > *i > ji > *2 , #2 e {1 , . . . , n} , ii 7^ ji , ii ^ h ■ 
Under Markovian assumption we have much stronger result: 

sP(t,n) = ——, \<k<n + l, 
fc n + 1 

*") = (n+lHn + 2) ' + 

*S M = (n+l)(n + 2) ' i,ie{1 '-' n + 1} ' 
Note that 2q 2 = bq and 

£ 2 (^) U Gn(n) = {G E Q(n) : n x > 0, n 2 > 0, n x + n 2 < n} . 

Hence we just proved the following statement. 
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Lemma 4. In the Markovian case 

E(V(x(t)\U t = n) = K n N V^ + ^-[J2K n N 2 )^+ (18) 



+ 



n + 

t 



\n 2 =0 



n + 



T E K N k N~ n2 rf( ° )( l0 + (19) 



vn2=0 
t 2 



(n + 1 (n + 2 



E n 2 ^ 1 H- (20) 



\,ni+7i2<n 

Now we are going to evaluate EV(x(t) by averaging on Yl t : 
EV(x(t) = E(E(V(x(t)|n t )) = 

= E^r^E^wi 11 ^^)- ( 21 ) 

n=0 

To do this we need two technical lemmas. 

Lemma 5. For any finite- dimensional matrix A the following identities hold 

E# = e -*Nt(id-A) ( 22 ) 
t A Ut+1 1 

E-Jr = T- (e~ 5Nt(Id - A) - e-^Id) (23) 

+ 1 Oat 

E t2AUt+2 = 1 ( e -8*W-A) Id+5 N tA \ 

t (n, + i)(n t + 2) 5% \ e ) l/4j 

where (U t , t > 0) is the Poissonian process of intensity S^. 

Proof of Lemma is very easy. We omit it. 
Lemma 6. Leta\,a2 G (0,1). Consider 

t/i(a l5 a 2 ) := E ^afa^ 2 , 

* n 2 =0 

£Ma 1; a 2 ) := E — ^ — £ afa- 

rai+n 2 <n t 

where (H t , t >0) is a Poissonian process of intensity 5^ > 0. Then 
for ai ^ a 2 

U^ai, a 2 ) = Sjf 1 (ai - a 2 y l (e^ 1 " 01 ^ 1 _ e -(i-a 2 )%^ 

p — (1— ai)<5jv* p — (1— a2)<5jv* 



/ r » .. .. \ 1 1 — «1 1 — a 2 



5 N ■ U 2 {a 1 ,a 2 ) 



(1 - ax)(l - a 2 ) ai - a 2 
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and for a\ = a 2 = a 



U\(a, a) 
5 2 N ■ U 2 (a,a) 



te 



-(l—a)5pft 



(1-a) 2 \(l-a) 



5 Art 



1 - a 



-V \ e -(l-o)M 



Proof of Lemma El will be given in Section EJ 

Now we proceed with evaluation of (}2"Tj) . Putting A = in d23) and applying 
averaging (1211) to the first summand in (|18l) we get 



E (K"*V (0) ) =e -W(H-Jfff) E y(o) = e ^JVEF (0) , 

since Kn = L/v (see notation ( Tl2|) -( 1T3|) at the beginning of Subsection l4.2p . 

Similarly, using (T23l) for A = Kn and A = Id , we obtain 



v n 2 =0 



E n^i(5>"h = E n^l( Id -^ _1 (id-^ +1 ) qi 

= (-Z^)" 1 (Id - e^') qi . 



To find expectation of summands ( [19]) and (120]) we shall analyze spectrum 
of the operator L^v : IR 2 — )■ M 2 and then apply Lemma El It is easy to check 
that operator Ln has two different eigenvalues X^n and A 2i jv, corresponding to 
eigenvectors e.\ jy and e 2 ^, 



Ai. 



iV 



-ajv, A 2j at = -2 (ojat + 



(25) 



1 

2a N 



\ I 1 

2Ai ) .w 



\ a N + 2(3 N ) 



(26) 



\ A 2; jv — Ai^ / 

where oln and are the same as in ( Fl2|) . Hence actions of the operators k^ld 
and of the vector e^jy are identical, 

/cArIdei i7V = K N e hN = ae hN , a = 1 - ^ajv = 1 + 

but their actions on the vector e 2 jv are different: 

k N Ide 2tN = aie 2 ,AT, 

K N e 2tN = a 2 e 2tN , a 2 = 1 - 25^ (a N + f3 N ) = 1 + S^X^n , a x ^a 2 . 



0-1 = 1- £jv «iv = 1 + Sj/X^n , 
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In the basis euv, e 2) N we have 
2a n 



1 


2Ai. 



ei.Ar + w 2 ,Ne 2 ,N , where 



W 2 ,7V = = -T -T , (27) 

Substituting this decomposition into (|T9|) - (l20|) . we apply Lemma |6] to calculate 
expectations of (TT9j) and (120]) separately on each linear subspace (e^jv) and (e2,jv). 
Finally, noting that (1 — dj) <5at = —X^n we come to the following statement. 

Proposition 2. T7ie functions -Rjv(t) and D^it) can be given in the following 
explicit form: 



R N {t) 

D N (t) 



Rn(0) 
D N (0) 



-LnY 1 (ld-e LNt ) 



a" 
2a 2 



+ (v-r)d N (0) te Xl < Nt e hN + 



^2,N — Ai t 



W 2 ,N e 2 ,N + 



N 



+ 



+ 



-,Ai wt 



X 2 



A 2 jv ^1,JV 



,f-r) eijv + 



AitvA 



+ 



A 



2,7V 



3^1, JV* ' 

Ai,jv 



2,/V 



N 



\ 



w 2>N (v -r) 2 e 2:N , 



J 



where d^(t) := Ed(x(t)) with d(-) defined in 

Remark 5. Note that many parameters in the above formula depend on N. 
in the limit N — > oo these dependencies become rather simple: 



T 1 M 1 ( ~ a ° 

Ln ~N M ~ n{ 2a -2(a + /3) 



1 

2a 



e 2 N 



N 



a 

~N 



2(a + /3) 

A2,iV ~ 77- ) W 2,N 



a + 2/3 
2a 



But 





1 ' 



N 1 a + 2/3 

Remark 6. Proposition [2] can be rewritten in the following way 



i«jv(t) 



(7 

2a 2 



+ 



+[v - r) t djv(O) flfi (Ai jA rt) ei,jv H - ; — w 2,n e 2 ,7v 



(A 2) at — Xi,n) t 

+(v-r) t ( g 2 (\i,Nt) e 1>N H — w %N e 2 , N 

[A 2 ,N — M,N) t 



+ 
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e y - 1 

where gAy) = e v , g^iy) = • Using (1270 we easily transform this represen- 

y 

tation as follows 



2 

la 2 



+ (v-r)td N (0) U(X hN t) f J J + M9i,t) f J )) + 
+ {v-rft 2 (g' 2 {\ 1>N t) ( J ) + F N (92,t) ( J ) 



where 

F N {g,t) = — 77 : — — g (Xi, N t) w 2 ,n ■ 

\ [M,N — M,n) t J 

For the function d^(t) we also have explicit formula. 

Proposition 3. The function d,N{t) has the following form 

d N (t) = d N (0)ex V (-£t) + (l -exp (V ~* )N (28) 

and hence is a solution to the following equation 
d ex 

-d N (t) = --d N (t) + (v-r). (29) 

Proof of this proposition can be obtained by the same method as the proof 
of Proposition [2] but the corresponding reasonings are much more shorter and 
simpler. So we omit the proof of Proposition [3] and simply refer to [T8] where 
similar statement was presented in full details. 

The following statement follows from Proposition [2] by direct calculation. 

Corollary 1. (iijv(t), Djsf(t)) is a solution of the following system of differential 
equations: 



d_ f R N (t) \ / -a/N \ / R N (t) 

dt 1 D N {t) J I 2a/iV -2(a/N + P/(N-l)) J \ D N {t) 



+ 



2(u - r)dtf(t) \ / a 



r 2 



* o y^V 2 ^ 2 

where the function d]y{t) is the same as in [2$) and l[29]) . 
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4.5 Study of the asymptotical behavior 

In this subsection we obtain all results on asymptotic behavior which were pre- 
sented in Section [21 

Proof of Theorem Q] 

Here iV is fixed and t tends to infinity. Since matrix has two different negative 
eigenvalues one easily concludes that in explicit representation of Proposition |2] all 
terms containing e Livi or e Xi ' Nt go to zero as t — > oo. So the limit is equal to 

R N (oo) \ , r N _i / o 2 \ . , 2 / 1 1 
D N (oo) ) = { ~ Ln) {2^) + ^ {^ ei > N+ \^ W2 > Ne2 ' N 

Matrix (— L^) 1 can be calculated explicitly: 

(_T \-l _ f "JV 1 

Using ( |26|) and ( l2Tj) after some algebraic transformation we get 



( — \ 2 2 ( — \ 2 

R N (t) R N (oo) = a 2 a~ N l + l \ 2 T) = — N + {V P N 

A| N a a 



r2 



j „ . , 2 a 2 2(v-rf 2 a 2 A7 (u - r) 2 Ar2 
a at + Pn M,nM,n a + p a ■ (a + p) 

This proves Theorem [TJ 
Proofs of Theorems [2] and [31 

At the beginning we study the first time scale, namely, we assume that 

t N ->• 00 but t N /N ^ as iV — )• 00 . (30) 
To do this we shall use the next lemma. 

Lemma 7. Consider a function g — g(y) which assumed to be analytical in some 
neighborhood of the point y = 0. Assume also that g"(0) 7^ 0. Then for any 
sequences {^/tv} and {z^} tending to as N — >■ 00 we have 

, f n (zjv) - 9 (2/Jv) 

P (yjv) e 1JV H W2,Ne 2 ,N = 

zn - Vn 

= 9 (vn) ( q J + Q /(0) (z N - y N ) w 2 ,n + (z N - Vn) O (\z n \ + \y N \)^ y J 
where the notation = 0{b^) means, as usual, that limsup 

N 



< +OO. 
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Proof of this lemma is a straightforward calculation. 

Since the both A^jv and A2.7V are negative and condition (J7J) holds we remark 
that the first summand in the r.h.s. of Proposition [2] is bounded. Under assump- 
tion ([30]) 

Ai,7vtjv — > and X 2j Nt N — > . 

So the second summand has a rather simple behavior: 



2 \ / 2 



W < Id - eL "'") { l> ) ~ { * + ( W ) ( ~ • < 31 > 

The next summands of the representation in Proposition [2] demands very careful 
analysis. We use Remark [6] and Lemma [7] for y N = A^jytjv anc ^ ^iv = ^n^n- 

^i(Ai,Arrj ei jA r H 77 r r- W2,tv e 2 .7v = 

= em ( ) + Q ^ ~ ^ ^ + ^ zn ~ y ^ ° ( tN / N> )) ( 1 

= (1 + 0(X ljN t N ) ) ( J ) + ( |Ai,iv| t N + O (t%/N 2 )) f J 

l + 0(t N /N) \ 
\\i, N \tN + 0(t 2 N /N 2 ) J [6Z) 

where we have use identity (z N — y N ) w 2 ,n = 2 (Ai^l t N (see and ( 127)) ). Sim- 
ilarly, 

//\ a £2 (A 2 ,jvO - fl-2 (Ai,^) 

{*2,N — M,N) t 

= 1 ~ (1 ~ 2 ^ N)eVN ( J ^+(l-l(z N -yN)w 2 , N +(z N -y N )0(t N /N)\ ( J 

i + 0(Ai^)^ ( J ) + Q|A liJV |tv + 0(4/iV 2 )^) ( J 

l + 0(t N /N) \ 
l\\i, N \tN + 0(t%/N 2 ) J W 

Putting (l3T|) - (l33|) into formula of Remark |6] we get 

R N (t N ) = aH N + 0(t N /N) + (v — r) d N (0) t N + (v — r)t N 0{tN /N) + 
+ 1 -(v-r)Hl f + (v-r)HlO(t N /N) 
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D N (t N ) = 2o 2 t N + 0(t N /N) + 

+(v - r) d N (0) \X 1;N \ t% + (v- r)t N (t 2 N /N 2 ) + 

+ - r) 2 |A 1)JV | t% + (v- r)H 2 N (t 2 N /N 2 ) 

In the case v — r = these formulae immediately imply item PI of the Theorem [2j 
Consider the case v — r ^ . Since £jv oo terms containing (v — r)tjv are 
asymptotically smaller than corresponding terms containing {y — r) 2 t 2 N . Moreover, 
\^\,n\ t 2 N = t N 0(t N /N). Hence 

R N (t N ) = a 2 t N + 0(t N /N)+ (v-r)d N (0)t N + 

+ l -{v-r) 2 t 2 N + {v-r)H 2 N 0{t N /N) 



D N (t N ) = 2a 2 t N + O(t N /N) + (v-r)d N (0)t N O(t N /N) + 

+ ~ r) 2 \\i, N \ t% +(v- r)H 2 N {t 2 N /N 2 ) . (34) 

So our next conclusion is that under assumptions fl30l) and v — r ^ 
RnM ~ \{v-r)H 2 N . 

To analyze Dn^n) we should compare and IAi^jvI ^jv- Recalling that |Ai,jv| = 
a 7v = oc/N we split the time scale (1301) into three subscales 

— Y U, — p= — > Ci, ci > (J, and — = — > oo but — — >■ 0. 



ViV ViV viV N 

Considering (154]) on each subscale we get that 

• if 0, then ~ 2<7 2 t JV , 

• if -^L ->■ ci, ci > 0, then 

D N (t N ) ~ 2a 2 t,v + g«(v - r) 2 4/iV 



~ 2a 2 ciVN + -a(u - r) 2 c?ViV ~ ^2a 2 + -a(u - r) 2 c 2 ) t/v 
if ^ oo but ^ ->■ 0, then D N (t N ) ~ |a(u - rf t%/N. 
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Next let us prove Theorems [2] and [3] for the time scale 

t N /N -> c> 0, N oo. (35) 
Taking into account Remark |5] and assumption (135]) we get 

£atLtv - > cM, Xi,NtN —&c, A2,Artiv — > —2(a + j3)c 

where M is defined in (jHJ). If is an analytic function then (/(A^jvt/v) = g'(-otc) + 
o(l) and Fjv(5', i/v) = # (#, c) + o(l), where 

fg(-2(a + P)c) -g(-ac) \ ( 2a 

H (9,c) := — — rr — ■ (-ac) 



-2(a + l3)c + ac v '/ \ a + 2/3 

Using Remark we have 

RN(tN) \ ( „i/r\-\ (t A „cAf\ ( ° 2 



sX)j^ ( - cMri(id - ecA n^j^ + ^ )+ 

+ (v-r)t N Ed N (0) Ui-ac) ( J ) + H(g u c) (J V 
+ (^-r) 2 4f(^(-ac) + o(l)) ( M + c) + (l)) / ° 



o(l) 
o(l) 



+ 



Now if u = r we easily get item P2 of Theorem |5J In the case u/rwe obtain 
item P2 of Theorem [3] with 

h R {c) = g' 2 (-ac), 
h D (c) = H(g 2 ,c), 

where g% is the same as in Remark [6j Now we are able to get explicit forms of the 
functions hf>(c) and /id(c): 

j , \ 1 - (1 + ac)e- ac 
Mc) = 3^ ' 



crc z 

_ 2 1 — (1 + ac)e" 



a (a + 2pJ 

2c~ 2 a / 1 1 /e" ac e - 2(a+/3)c 



a + 2/3 \2(a + /3)a a + 2/3 \ a 2(a + /3) 

The study of the time scale tw/N — > oo is very similar to the proof of Theorem [1] 
(see above). So we omit details. 
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4.6 Time scale analysis 

In this subsection we derive some corollaries from our main theorems. We examine 
behavior of the functions Ji^^jv) an d -Djv(^iv) on special time scales — sN 1 
where 7 > and s > 0. Summing up Theorems [TH3] and Subsection 14.51 we get 
the following general statement: Assume that iV — > 00 and assumption ((7|) holds. 
Then 



C D {s,i) 



D N (sNi) 

where all functions Cr, Cp, ip R , ipD are positive. Moreover, these functions can be 
calculated explicitely: 



Case 1: v = r — zero skew. 

^(7) = min(7, 1), C R (s,i) 

^d(t) = min(7, 1), C D (s, 7) 
where 




Ir(s) 
Id(s) 



g 2 {-as) , 

a g 2 {-as) + 2(5 g 2 (-2(a + 0)s) 



a + 2/3 

and the function g 2 is the same as Remark [6j 
Case 2: v 7^ r — nonzero skew. 



\{v — r) 2 s 2 , 
(v - rfs 2 h R (s) 
(v — r) 2 



^r{i) = min(27,2), 



2' 



7< 
|<7<1, 
7> 1, 



a' 



2a 2 s, 



7 < 1 
7 = 1, 

7 > 1, 



2a 2 s + \a{v 



r) 2 s 3 , 



7< h 
1 = \ 



2 „3 



|<7<1, 



where h R and ho are the same as in (136|) . 



— r) 

(u -r) 2 s 2 h D (s), 
k (u -r) 2 «- 1 (a + 7>1. 



2 

7 = 1 
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How to interpret and to explain the above results? The choice of the time scale 
tjy = sN 1 means that we consider a new time unit which is equal to iV 7 units of 
the physical time t. Then paramemer s is the time measured in new units. We 
see that behavior of the stochastic system x(t) = (x\(t), . . . , XN+i(t)) for large 
N is different for different values of 7. We recall (Section EJ) that R^lsN 1 ) and 
Dn{sN j ) are squares of clock synchronization errors. The above corollaries show 
existence of several phases in the evolution of the network. In the case v = r there 
are three such phases, they were discussed in Remark EJ 

In the case u^ron "small times" (7 < 1) the function Rn(sN' j ) is propor- 
tional to the square of t^. We can imagine that on these times effect of clock 
synchronization is negligible and the asymptotical value R^isN 1 ) depends only 
on summands of the form 

(xx(t) - Xj (t)) 2 = (xx(0) - Xj (0) + (r - v)tf ~ (v - rft 2 . (37) 

For the scale 7=1 the total number of synchronization messages from the server 1 
is of order N, their influence become important and we observe a slowdown: for 
s — > 00 

s 2 h R (s) | < +oo- 

On the scales 7 > 1 we see the result of full strength competition between syn- 
chronization jumps and the desynchronization generated by motions of the client 
nodes 2, . . . , N + 1. This is the final synchronization phase. 

In in the same case v 7^ r the behavior of the function D^lsN 1 ) is much 
more interesting. It appears that the "small times" (7 < 1) are splitted into 
three subphases: 0<7<|,7=| and | < 7 < 1. Recall that the function 
describes the internal inconsistency of the client's clocks X2, ■ ■ ■ , a^jv+i- On the first 
subphase (0 < 7 < |) influence of the server node 1 is negligible in comparison with 
a "noise" produced by identical client nodes, hence we see the same asymptotics 
2<t 2 s iV 7 as in the Case 1 (see also asymptotics in [18j). On the scale 7 = | these 
two forces are of the same order N 1 / 2 and Cd{s, |) = 2a 2 s + |a(f — r) 2 s 3 . On the 
third subphase (| < 7 < 1) the server node 1 dominates over the free dynamics 
of the clients. Surprisingly, that for scales ~ < 7 < 1 the influence of the server of 
the accurate time produces a desynchronization (but not synchronization) of the 
client's clocks. The explanations is the following one: on these scales only small 
part of client nodes had interacted with the node 1 till the time sN 1 . The most 
part of client nodes (of order N) still have no idea about the server's clock x\. 
Adjusting of a client's clock to the value x± brings a large number of summands 
of the form (157)) to the function D^. On the scale 7 = 1 the number of nodes 
interacted with the server 1 is of order N and again we observe a slowdown of 
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desynchronization. The time scales 7 > 1 correspond to the final synchronization 
phase. 



5 Appendix 

Proof of Lemma [6], Let a\ 7^ It is straightforward to check that 

n 

a?ar n2 = (ax - a^ 1 « +1 - a^ 1 ) . 

n 2 =0 

Using f l2~5j) for A — a,\ and A = a% we find U\{a,\, a<i). 

Now let us calculate Ui{a\,a-i). After some simple algebra we have identity 

_1 1 a? +2 a?+ 2 



,, 1 „.j 1 — ai 1 — 1 — ai 1 — 02 



ai — a 2 ai — a 2 

ni > 0, n2 > 

ni + ri2 < n 



Applying ([21]) for A = 1, A = ai and A = a 2 , we get that ^(ai, a 2) multiplied by 
5 2 N is equal to 

1 1 



1 - ax l-a 2 ( l _ 1 + <M 



ai — a 2 



1 — ai 1 — a 2 



ai — 02 

By direct transformations and cancellation of terms this form can be reduced to 
the following one 

e -<5jvt(l-ai) _ j g -(5jvt(l-a2) _ I e -8 N t(l-ai) e ~8 N t(l-a2) 



1 — ai 1 — «2 1 1 — ai 1 — a2 



ai — a 2 (1 — ai)(l — a 2 ) ai — a 2 

This proves the statement of Lemma [6] for a% 7^ a 2 . 

The case ai = a 2 is simpler than just considered case a% 7^ 02- So we omit 
details here. Note, that explicit expressions for Ui(a,a), i = 1,2, correspond to 
formal limits of £/j(ai, 02) as ai — > a, 02 — > a. □ 
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6 Conclusion and future work 

We proposed a basic probabilistic model of clock synchronization in large WSNs 
interacting with an accurate time server. It was shown that in large networks 
(N — > oo) there exists several time scales t = t jy of qualitatively different collective 
behavior of the network. In other words, a large network passes different phases 
on its road to synchronization. The phase of effective synchronization is the most 
interesting among them. For our basic network we give a detailed description of 
this phase. Moreover, explicit formulae obtained in Theorems [TH3] provide keys 
to future analytical study of various optimization and performance evaluation 
problems related to WSNs. 

We believe that obtained results about existence of several different phases in 
evolution should take place for more general classes of large networks. Mathemati- 
cal tools used in the present study will work also for general nonhomogeneous time 
synchronization model. Of course, in that case we cannot expect such short and 
explicit results as in Section [3j Successful examples of synchronization studies for 
some special nonhomogeneous (or weakly nonhomogeneous) systems can be found 

in [sinn Ea Eg. 

Our basic model can be developped also for more general and realistic assump- 
tions about message sending algorithms. We plan to get rid of condition ([3]) to be 
able to consider arbitrary distributed intervals between messages. The theory of 
general random flows [2J can be useful here but, unfortunately, the corresponding 
probabilistic model will be non-Mar kovian (similarly to [19]). The white noise 
assumption in the clock model (jTJ) is confusing for the following reason: it con- 
tradicts to a general clock modelling principle that time never run backward (see 
page 313 in [30j). To justify our choice we note, first, that from the mathematical 
viewpoint this assumption is convenient but not necessary and it will be removed 
in future papers, and, in the second place, this is an usual assumption for many 
modern clock models (see [HI [321 E]) 

We hope that presented approach combined with other methods will be useful 
for analysis of many practical distributed systems. 
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